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QCD predictions for moments of multiplicity distributions are compared with 
experimental data on e^e~ collisions and their two-NBD fits. Moments of the mul- 



I tiplicity distribution in a two-NBD model for 1.8 TeV pp-collisions are considered. 

Three-NBD model predictions and fits for pp at LHC energies are also discussed. 
. Analytic expressions for moments of hybrid NBD are derived and used to get insight 

into jet parameters and multicomponent structure of the processes. Interpretation 
^ ' of observed correlations is proposed. 

m ' 

Multiplicity distributions are the integral characteristics of multiparticle production 
processes. They can be described either in terms of probabilities Pn{E) to create n 
^ ' particles at energy E or by the moments of these distributions. It has been found that 

O ' their shapes possess some common features in all reactions studied. At comparatively 

low energies below tens of GeV, these distributions are relatively narrow and have sub- 
Ph! Poissonian shapes. At energies about 20 GeV for e~^e~- annihilation and 30 GeV for 

f^l pp (and pp) interactions, they can be well fitted by the Poisson distribution^. At higher 

energies, the shapes become super-Poissonian, i.e. their widths are larger than for Poisson 
distribution. They increase with energy and, moreover, some shoulder-like substructures 
appear. 

Their origin is usually ascribed to multicomponent contents of the process. In QCD 
5^ I description of e"*"e~-processes these could be subjets formed inside quark and gluon jets 

(for the reviews see, e.g., [1, 2]). In phenomenological approaches, the multiplicity distri- 
bution in a single subjet is sometimes approximated by the negative binomial distribution 
(NBD) first proposed for hadronic reactions in [3]. For hadron-initiated processes, these 
peculiarities are also explained by the multicomponent structure of the process. This is 
either multiladder exchange in the dual parton model [4, 5] , varying number of clans [6] or 
multiparton interactions [7, 8]. These subprocesses are related to the matter state during 
the collision (e.g., there are speculations about nonhomogeneous matter distribution in 
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^This difference of thresholds can be attributed to lower energy spent in hadron collisions for new 
particles creation (due to the so-called leading particle effect) compared to e^e^- annihilation. 
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impact parameters [9], not to speak of quark-gluon plasma [10] behaving as a liquid [11] 
etc). 

Such evolution of the multiplicity distributions can be quantitatively described by the 
energy behaviour of their moments. These moments reveal the correlations inherent for 
the matter state formed during the collision. Similarly to virial coefficients in statistical 
physics they can tell us about the equation of state of this matter. To introduce them, 
let us write the generating function of the multiplicity distribution as 

oo 

G{E,z) = Y.Pn{E){l + zr. (1) 

n=0 

In what follows, we will use the so-called unnormalized factorial !Fq and cumulant ICq 
moments defined according to the formulae 
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They correspondingly define the total and genuine correlations among the particles pro- 
duced (for more details see [12, 2]). These cumulant moments could be considered as the 
direct analogies of virial coefficients of statistical physics since both arc related to genuine 
(irreducible) correlations. In particular, the first moments describe the mean multiplicity 

in) 

^i = /Ci = (n), (4) 
and the second moments are related to the dispersion D of the distribution P„: 

/C2 = ^2-(n)' = i^'-(n). (5) 

The higher rank moments reveal other asymmetries of distributions such as skewness etc. 
Since both and fCg strongly increase with their rank and energy, their ratio 

Hg = ICg/Tg (6) 

first introduced in [13] is especially useful due to partial cancellation of these dependences. 
The factorial moments TgS are always positive by definition (Eq. (2)) while the cumulant 
moments /C^'s can change sign. Again, let us recall that the changing sign second virial 
coefficient in statistical physics implies the liquid state with the Van der Waals equation 
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corresponding to repulsion at small distances and attraction at large distances. Cooper 
pair formation is also related to similar behaviour of correlations. 

Here, we compare QCD and NBD approaches to the description of multiphcity distri- 
butions. We argue that Hg values are more sensitive to minute details of the distributions 
than their direct chi-square fits and reveal differences between proposed fits of e'^e~ and 
PP {pp) processes. Some estimates for LHC energies will be provided. 

The generating functions for quark and gluon jets satisfy definite equations in pertur- 
bative QCD (see [14, 2]). It has been analytically predicted in gluodynamics [13] that 
at asymptotically high energies the Hg moments are positive and decrease as but 
at present energies they become negative at some values of q and reveal the negative 
minimum at 

qmin^T + 0.5 + 0(70), (7) 

/ii7o 

where hi = b/SN^ = 11/24, b = llA^c/3 - 2n//3, 7^ = 2Ncas/n, as is a coupling 
strength, N^, nj are the numbers of colours and flavours. At 2"° energy ~ 0.12, 
and this minimum is at about g ^ 5. It moves to higher ranks with energy increase 
because the coupling strength decreases. Some hints to possible oscillations of Hg vs q at 
higher ranks at LEP energies were obtained in [13]. Then the approximate solution of the 
gluodynamics equation for the generating function [15] agreed with this and predicted the 
oscillating behaviour at higher ranks. These oscillations were conflrmed by experimental 
data for e+e" and hadron-initiatcd processes flrst in [16], later in [17] and most recently 
in [18]. The same conclusions were obtained from exact solution of equations for quark 
and gluon jets in the framework of fixed coupling QCD [19]. The physics interpretation 
of these oscillations as originating from multisubjet structure of the process is related to 
the (multi)fractal behaviour of factorial moments, found also in QCD [20, 21, 22]. The 
asymptotic disappearance of oscillations can be ascribed to the extremely large number 
of subjets at very high energies. 

A recent exact numerical solution of the gluodynamics equation in a wide energy 
interval [23] coincides with the qualitative features of multiplicity distributions described 
above. In terms of moments they correspond to the values of Hg changing sign at each 
subsequent q (with H2 < 0) at low energies (narrow shapes^), the approach of Hg to 
zero at the Poisson transition point about 20 GeV for e'^e~ processes, and the positive 
second moment H2 with oscillations of higher rank cumulants at Z*^ which disappear 
asymptotically. At Z°, the first minimum appears at ~ 5. This confirms earlier exact 
QCD results [21] at Z". It mows to higher ranks with a steadily decreasing amplitude 

^Narrow distributions always have such cumulants as shown, e.g., in [2]. 
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when energy increases. The only free parameter is the QCD cut-off, which is however 
approximately fixed by the coupling strength and does not strongly influence the results. 

In parallel, the NBD-fits of multiphcity distributions were attempted [6, 25]. The 
single NBD-parameterization is 

where F denotes the gamma-function. This distribution has two adjustable parameters 
ni{E) and ki{E) which depend on energy. Such a formula happened to describe low 
energy data with negative values of ki that corresponds to binomial fits. At the Poisson 
transition point = 0. The parameter ki becomes positive at higher energies. However 
the simple fit by the formula (8) is valid till the shoulders appear. In that case, this formula 
is replaced by the hybrid NBD which combines two or more expressions like (8). Each of 
them has its own energy dependent parameters n^, /cj. These distributions are weighted 
with the energy dependent probability factors ctj which sum up to 1. Correspondingly, 
the number of adjustable parameters drastically increases. 

A single NBD (8) has positive cumulants for ki > (ICq = r{q)nl/ki^^) and thus the 
positive Hq = r(g)r(A;i + l)/r{ki + q). For hybrid NBD, the negative Hg can exist. The 
traditional procedure to calculate higher rank moments is by the iterative relations 

The strong compensations are inherent in Eq. (9). This calls for high accuracy of numer- 
ical calculations. More important, the formula does not give any direct insight into the 
physical reasons for such compensations. Therefore, it is instructive to write the analytic 
formulae for moments of hybrid NBD which provide clear interpretation of negative val- 
ues of cumulants. We have derived these expressions for the two-NBD parameterization 
(2NBD) given by a sum of two expressions like (8) with two sets of adjustable parameters 
ni, /ci, n2, k2 weighted with energy dependent factors a and 1 — a correspondingly. 2NBD 
describes the process with two independent NBD-components of mean multiplicities rii 
and widths ki created with probabilities a and 1 — a. The factorial moments for any rank 
q are given by the simple formula 



r{k,) ky^ V{k2) kl 



•Fq^OL ^ -^\{\-a) ^ (0 < a < 1). (10) 
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The cumulant moments are more complicated and should be calculated separately for 
each rank. The first 5 moments are 




„ an? (1 — a)nl ,^ , , ,0 

Ki K2 




/C3 = ^ + J'^''' + a{l - a){n, - n^Mnl/k, - nl/k^) + (1 - 2a){n, - n^fl (13) 



/C4 = ^ + f""^ + a(l - a)[(ni - n^fil - 6a(l - a)) + Il(n2/A;i - n^/A;2) - 



8nin2(ni/A;i - 712^2)^ + 6(1 - 2a){ni - n2f{n\/ki - n^/A;2)], (14) 




For q; = or 1 they reduce to one-NBD formulae with one of the first two terms surviv- 
ing. It is always positive for positive k^. Therefore, as expected, individually considered, 
the distributions show no oscillations. For 2NBD, there is a symmetry in replacing in- 
dices 1 to 2 together with a to 1 — a. Negative /C2 can be obtained only ii ki < 0. For 
positive ki one always gets positive /C2. Its value depends on the difference rii — 7^2. /C3 
can become negative depending on the values of the last two terms. The cancellations in 
the expressions (13)-(15) are not so drastic as in Eq. (9), especially for large q, because 
the leading contributions to Hq are strongly decreasing with q there and not of the order 
of 1 as in (9). Therefore they do not require very high precision and, moreover, clearly 
display the origin of each term and its dependence on fitted parameters. 

Actually, five moments determine quite well the shape of the distribution if they are 
calculated with high enough accuracy. Since these shapes arc qualitatively similar in 
different reactions, it is especially instructive to compare their Hq moments. In Table 1 
the Hq moments for e~^e~ annihilation at are shown. Their values according to the 
solution of the gluo dynamics equations [23] are in the first column. In the second and 
third columns, the experimental results of L3 collaboration [18] are represented for full 
phase space (FPS), correspondingly, with all measured multiplicities included and with 
some very high multiplicities truncated (because of large error bars). Next follow Hq 
values restored from 2NBD fits of OPAL and DELPHI results done in [26]. 
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Table 1. 





QCD 


L3, untr 


L3, tr 


0PAL,2NBD 


DELPHI,2NBD 




3.9E-2 


(4.42±0.11)E-2 


(4.41±0.10)E-2 


4.4E-2 


3.1E-2 


Hs 


7.4E-3 


(7.40±0.38)E-3 


(7.20±0.35)E-3 


7.4E-3 


2.6E-3 


H4 


4.0E-4 


(9.69±2.56)E-4 


(7.17±1.42)E-4 


4.9E-4 


7.4E-5 


H5 


-2.2E-4 


-(1.30±1.59)E-4 


-(3.95±0.53)E-4 


-2.4E-4 


-7.3E-5 



The overall agreement is rather good. Quite impressive is the fact that in all cases the 
fifth cumulant moment is negative. However, somewhat surprising is the difference of the 
theoretical and experimental widths {H2 values). The widths are determined quite pre- 
cisely both experimentally and theoretically. The only reason to which such disagreement 
could be ascribed is the incomplete treatment with quarks omitted in [23]. The more 
complete approach will shed more light on this problem. The even stronger disagreement 
with DELPHI data is probably related to the special selection of events there. Further 
analysis is needed. 

The comparison of e+e^ and pp {pp) data turns out especially interesting. While both 
show qualitative similarity of the shapes of multiplicity distributions, the corresponding 
Hq values are quite distinctive. In Table 2 we show Hg values for pp at 1.8 TeV (Tevatron) 
and interpolations to 14 TeV (LHC). The 2NBD fit at 1.8 TeV corresponds to the follow- 
ing parameters: a = 0.62, Ui = 30, n2 = 61.6, ki = k2 = 7, which are approximately 
equal to average values for 2. A model considered in [6]. However, even this extreme 
model underestimates high multiplicities and, therefore, Hg values in the Table should be 
treated as lower bounds to experimental ones, which are unknown, unfortunately. The 
extrapolated values at 14 TeV have been calculated using the parameters of 3NBD fits 
and Pythia model considered in [25]. 



Table 2. 





2NBD fit, 1.8TeV 


3NBD fit, 14TeV 


Pythia, 14TeV 


H2 


0.2279 


0.8754 


0.4224 




0.0988 


0.9703 


0.3387 




0.0414 


0.9737 


0.2683 




0.0120 


0.9712 


0.1877 



Quite impressive are much larger values of Hg in hadron-initiated reactions (Table 2) 
as compared to e+e~ results (Table 1). They strongly increase with energy. Moreover, 
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the drastic difference is clearly displayed by Hq between 3NBD interpolations and Pythia 
at 14 TeV. This demonstrates extremely high sensitivity of Hq analysis because both 
approaches provide the similar two-shoulder structure of multiplicity distributions as seen 
in Fig. 2 of [6]. At 14 TeV, the predictions are given for full phase space. For the rapidity 
interval [t^I < 0.9 the Hq values become larger than those in Table 2. Hq for the 3NBD- 
model of [6] become almost indistinguishable from 1 (above 0.99). Pythia values increase 
about 1.4 times. No oscillations are seen at these high energies while they are present at 
energies below 1 TeV [16]. Surely, LHC experiments will give their decisive conclusion. 

To conclude, we have shown that Hq moments of the multiplicity distribution are 
extremely sensitive to minute details of its shape. They can resolve the differences between 
various fits even if those are not clearly seen in the traditional representation. For e+e~, 
slight disagreement on theoretical and experimental widths is embarassing and must be 
further studied. For hadron- and nuclei-initiated reactions, Hg values are much larger than 
in e~^e~ . RHIC and LHC data are awaited for better insight. The energy dependence of 
Hq and of the relative weights of various NBD components can provide some hints on the 
matter state during the collision and its energy evolution. 

I am grateful to A. Giovannini and W. Metzger for correspondence. This work has been 
supported in part by the RFBR grants N 02-02-16779, 03-02-16134, NSH-1936.2003.2. 
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